{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Uncertainty calculation for model: WK2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# import the libraries\n",
    "import ee\n",
    "import pandas as pd\n",
    "import os\n",
    "import numpy as np\n",
    "import random\n",
    "from random import sample\n",
    "import itertools \n",
    "import geopandas as gpd\n",
    "from sklearn.metrics import r2_score\n",
    "from termcolor import colored # this is allocate colour and fonts type for the print title and text\n",
    "from IPython.display import display, HTML"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "#set the working directory of local drive for Grid search result table loading\n",
    "# os.getcwd()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# initialize the earth engine API\n",
    "ee.Initialize()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1 Load the required composites, images and settings"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "#definet the color pallette\n",
    "vibgYOR = ['330044', '220066', '1133cc', '33dd00', 'ffda21', 'ff6622', 'd10000']\n",
    "compositeImage =ee.Image(\"users/leonidmoore/ForestBiomass/20200915_Forest_Biomass_Predictors_Image\")\n",
    "compositeImageNew = ee.Image(\"projects/crowtherlab/Composite/CrowtherLab_Composite_30ArcSec\");\n",
    "unboundedGeo = ee.Geometry.Polygon([-180, 88, 0, 88, 180, 88, 180, -88, 0, -88, -180, -88], None, False)\n",
    "# generete the pixel area map\n",
    "pixelArea = ee.Image.pixelArea().divide(10000) # to ha unit\n",
    "# load the biome layer\n",
    "biomeLayer = compositeImage.select(\"WWF_Biome\")\n",
    "biomeMask = biomeLayer.mask(biomeLayer.neq(98)).mask(biomeLayer.neq(99)).gt(0)\n",
    "# load the mean maps for present and potential\n",
    "# load the carbon concentration map\n",
    "carbonConcentration = ee.Image(\"users/leonidmoore/ForestBiomass/Biome_level_Wood_Carbon_Conentration_Map\")\n",
    "# load the biomass density layers\n",
    "mergedAGB_PresentMean =  ee.Image(\"users/leonidmoore/ForestBiomass/WalkerMap/reprojected_Walker_map_1km\").unmask() \n",
    "mergedAGB_PotentialMean = ee.Image(\"users/nordmannmoore/ForestBiomass/RemoteSensingModel/EnsambleMaps/Predicted_WK2_Potential_density_Ensambled_Mean\").unmask()\n",
    "\n",
    "# define the standardized projection\n",
    "stdProj = mergedAGB_PresentMean.projection()\n",
    "# load the two forest cover layer for existing and potential forest\n",
    "presentForestCover = compositeImage.select('PresentTreeCover').unmask()# make sure it's in  0-1 scale\n",
    "potentialForestCover = ee.Image(\"users/leonidmoore/ForestBiomass/Bastin_et_al_2019_Potential_Forest_Cover_Adjusted\").unmask() # make sure it's in  0-1 scale\n",
    "\n",
    "# define the present and potential forest cover masks\n",
    "presentMask = presentForestCover.gt(0)\n",
    "potentialMask = potentialForestCover.gt(0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2 Calculate the present and potential AGB"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# check the difference of the two density maps\n",
    "potentialHigher = mergedAGB_PotentialMean.multiply(pixelArea).subtract(mergedAGB_PresentMean.multiply(pixelArea)).gte(0)\n",
    "potentialLower = mergedAGB_PotentialMean.multiply(pixelArea).subtract(mergedAGB_PresentMean.multiply(pixelArea)).lt(0)\n",
    "# replace the lower potential value by present biomass density value\n",
    "potentialAGB_Density = mergedAGB_PresentMean.multiply(potentialLower).add(mergedAGB_PotentialMean.multiply(potentialHigher))\n",
    "presentAGB_Density = mergedAGB_PresentMean\n",
    "# get the abs of present and potential AGB\n",
    "presentAGB_Abs = presentAGB_Density.multiply(pixelArea).multiply(presentMask).divide(1000000000)\n",
    "potentialAGB_Abs = potentialAGB_Density.multiply(pixelArea).multiply(potentialMask).divide(1000000000)\n",
    "\n",
    "# presentAGB_Abs_Sum = presentAGB_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                  geometry = unboundedGeo,\n",
    "#                                                  crs = 'EPSG:4326',\n",
    "#                                                  crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                  maxPixels = 1e9)\n",
    "# # print the estimation out\n",
    "# print(colored('The present AGB:', 'blue', attrs=['bold']),presentAGB_Abs_Sum.getInfo())\n",
    "# potentialAGB_Abs_Sum = potentialAGB_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                      geometry = unboundedGeo,\n",
    "#                                                      crs = 'EPSG:4326',\n",
    "#                                                      crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                      maxPixels = 1e9)\n",
    "# # print the estimation out\n",
    "# print(colored('The potential AGB:', 'blue', attrs=['bold']),potentialAGB_Abs_Sum.getInfo())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3 Calculate the Upper and Lower of present and potential AGB"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# load the present prediction lower and upper layer\n",
    "# this is the interval index, in %. divide by 2 and 100 to make it into the upper and lower boundary\n",
    "mergedPredictionInterval = ee.Image(\"users/leonidmoore/ForestBiomass/WalkerMap/reprojected_Walker_map_1km_Uncertainty\").unmask().divide(200)\n",
    "presentAGB_Lower = mergedAGB_PresentMean.subtract(mergedAGB_PresentMean.multiply(mergedPredictionInterval))\n",
    "presentAGB_Upper = mergedAGB_PresentMean.add(mergedAGB_PresentMean.multiply(mergedPredictionInterval))\n",
    "\n",
    "# get the upper and lower layer\n",
    "mergedAGB_PotentialLower = ee.Image(\"users/nordmannmoore/ForestBiomass/RemoteSensingModel/EnsambleMaps/Predicted_WK2_Potential_density_Ensambled_Percentile\").select(['lower']).unmask() \n",
    "mergedAGB_PotentialUpper = ee.Image(\"users/nordmannmoore/ForestBiomass/RemoteSensingModel/EnsambleMaps/Predicted_WK2_Potential_density_Ensambled_Percentile\").select(['upper']).unmask() \n",
    "\n",
    "mergedAGB_PotentialLower1 = presentAGB_Lower.multiply(potentialLower).add(mergedAGB_PotentialLower.multiply(potentialHigher))\n",
    "mergedAGB_PotentialUpper1 = presentAGB_Upper.multiply(potentialLower).add(mergedAGB_PotentialUpper.multiply(potentialHigher))\n",
    "\n",
    "# define the masks to mask the present and potential lower maps\n",
    "potentialAGB_Lower_Larger = mergedAGB_PotentialLower1.subtract(presentAGB_Lower).gte(0) # potential is larger than present mean\n",
    "potentialAGB_Lower_Smaller = mergedAGB_PotentialLower1.subtract(presentAGB_Lower).lt(0)\n",
    "# define the masks to mask the present and potential upper maps\n",
    "potentialAGB_Upper_Larger = mergedAGB_PotentialUpper1.subtract(presentAGB_Upper).gte(0) # potential is larger than present upper\n",
    "potentialAGB_Upper_Smaller = mergedAGB_PotentialUpper1.subtract(presentAGB_Upper).lt(0)\n",
    "\n",
    "# replace the lower potential value by present biomass density value\n",
    "potentialAGB_AdjLower = mergedAGB_PotentialLower.multiply(potentialAGB_Lower_Larger).add(presentAGB_Lower.multiply(potentialAGB_Lower_Smaller))\n",
    "potentialAGB_AdjUpper = mergedAGB_PotentialUpper.multiply(potentialAGB_Upper_Larger).add(presentAGB_Upper.multiply(potentialAGB_Upper_Smaller))\n",
    "\n",
    "# present lower and higher\n",
    "presentAGB_Lower_Abs = presentAGB_Lower.multiply(pixelArea).multiply(presentMask).divide(1000000000)\n",
    "presentAGB_Upper_Abs = presentAGB_Upper.multiply(pixelArea).multiply(presentMask).divide(1000000000)\n",
    "\n",
    "# abs potential lower and higher\n",
    "potentialAGB_Lower_Abs = potentialAGB_AdjLower.multiply(pixelArea).multiply(potentialMask).divide(1000000000)\n",
    "potentialAGB_Upper_Abs = potentialAGB_AdjUpper.multiply(pixelArea).multiply(potentialMask).divide(1000000000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # Calculate the present AGB lower\n",
    "# presentAGB_Lower_Sum = presentAGB_Lower_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                          geometry = unboundedGeo,\n",
    "#                                                          crs = 'EPSG:4326',\n",
    "#                                                          crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                          maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present AGB Lower:', 'blue', attrs=['bold']),presentAGB_Lower_Sum.getInfo())\n",
    "\n",
    "# # Calculate the present AGB upper\n",
    "# presentAGB_Upper_Sum = presentAGB_Upper_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                          geometry = unboundedGeo,\n",
    "#                                                          crs = 'EPSG:4326',\n",
    "#                                                          crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                          maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present AGB Upper:', 'blue', attrs=['bold']),presentAGB_Upper_Sum.getInfo())\n",
    "\n",
    "# potentialAGB_Lower_Sum = potentialAGB_Lower_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                              geometry = unboundedGeo,\n",
    "#                                                              crs = 'EPSG:4326',\n",
    "#                                                              crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                              maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential AGB Lower:', 'blue', attrs=['bold']),potentialAGB_Lower_Sum.getInfo())\n",
    "\n",
    "# potentialAGB_Upper_Sum = potentialAGB_Upper_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                              geometry = unboundedGeo,\n",
    "#                                                              crs = 'EPSG:4326',\n",
    "#                                                              crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                              maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential AGB Upper:', 'blue', attrs=['bold']),potentialAGB_Upper_Sum.getInfo())\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4 Calculate the Upper and Lower of present and potential Root and TGB"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# load the root shoot map\n",
    "rootShootRatio = ee.Image(\"users/leonidmoore/ForestBiomass/Root_shoot_ratio_Map\").unmask()\n",
    "rootShootRatioLower = ee.Image(\"users/leonidmoore/ForestBiomass/Root_shoot_ratio_percentile_Map\").select('lower').unmask()\n",
    "rootShootRatioUpper = ee.Image(\"users/leonidmoore/ForestBiomass/Root_shoot_ratio_percentile_Map\").select('upper').unmask()\n",
    "# \n",
    "presentRoot_Lower_Abs = presentAGB_Lower_Abs.multiply(rootShootRatioLower).mask(presentMask)\n",
    "presentRoot_Upper_Abs = presentAGB_Upper_Abs.multiply(rootShootRatioUpper).mask(presentMask)\n",
    "\n",
    "potentialRoot_Lower_Abs = potentialAGB_Lower_Abs.multiply(rootShootRatioLower).mask(potentialMask)\n",
    "potentialRoot_Upper_Abs = potentialAGB_Upper_Abs.multiply(rootShootRatioUpper).mask(potentialMask)\n",
    "\n",
    "presentRoot_Abs = presentAGB_Abs.multiply(rootShootRatio).mask(presentMask)\n",
    "potentialRoot_Abs = potentialAGB_Abs.multiply(rootShootRatio).mask(potentialMask)\n",
    "\n",
    "presentTGB_Abs = presentAGB_Abs.multiply(rootShootRatio).add(presentAGB_Abs)#.multiply(presentMask)\n",
    "potentialTGB_Abs = potentialAGB_Abs.multiply(rootShootRatio).add(potentialAGB_Abs)#.multiply(potentialMask)\n",
    "\n",
    "presentTGB  = presentAGB_Density.multiply(rootShootRatio.add(1))\n",
    "# density \n",
    "presentRoot = presentAGB_Density.multiply(rootShootRatio)\n",
    "presentRoot_Lower = presentAGB_Lower.multiply(rootShootRatioLower)\n",
    "presentRoot_Upper = presentAGB_Upper.multiply(rootShootRatioLower)\n",
    "\n",
    "potentialRoot_Lower = potentialAGB_AdjLower.multiply(rootShootRatioLower)\n",
    "potentialRoot_Upper = potentialAGB_AdjUpper.multiply(rootShootRatioLower)\n",
    "\n",
    "presentTGB_Lower = presentAGB_Lower.multiply(rootShootRatioLower).add(presentAGB_Lower)\n",
    "presentTGB_Upper = presentAGB_Upper.multiply(rootShootRatioLower).add(presentAGB_Upper)\n",
    "\n",
    "potentialTGB_Lower = potentialAGB_AdjLower.multiply(rootShootRatioLower).add(potentialAGB_AdjLower)\n",
    "potentialTGB_Upper = potentialAGB_AdjUpper.multiply(rootShootRatioLower).add(potentialAGB_AdjUpper)\n",
    "\n",
    "# presentTGB_Abs_Sum = presentTGB_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                  geometry = unboundedGeo,\n",
    "#                                                  crs = 'EPSG:4326',\n",
    "#                                                  crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                  maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present TGB:', 'blue', attrs=['bold']),presentTGB_Abs_Sum.getInfo())\n",
    "\n",
    "# potentialTGB_Abs_Sum = potentialTGB_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                      geometry = unboundedGeo,\n",
    "#                                                      crs = 'EPSG:4326',\n",
    "#                                                      crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                      maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential TGB:', 'blue', attrs=['bold']),potentialTGB_Abs_Sum.getInfo())\n",
    "\n",
    "# presentRoot_Abs_Sum = presentRoot_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                    geometry = unboundedGeo,\n",
    "#                                                    crs = 'EPSG:4326',\n",
    "#                                                    crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                    maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present Roots:', 'blue', attrs=['bold']),presentRoot_Abs_Sum.getInfo())\n",
    "\n",
    "# potentialRoot_Abs_Sum = potentialRoot_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential Roots:', 'blue', attrs=['bold']),potentialRoot_Abs_Sum.getInfo())\n",
    "\n",
    "# presentRoot_Lower_Sum = presentRoot_Lower_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present Roots Lower:', 'blue', attrs=['bold']),presentRoot_Lower_Sum.getInfo())\n",
    "\n",
    "# presentRoot_Upper_Sum = presentRoot_Upper_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present Roots Upper:', 'blue', attrs=['bold']),presentRoot_Upper_Sum.getInfo())\n",
    "\n",
    "# potentialRoot_Lower_Sum = potentialRoot_Lower_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential Roots Lower:', 'blue', attrs=['bold']),potentialRoot_Lower_Sum.getInfo())\n",
    "# potentialRoot_Upper_Sum = potentialRoot_Upper_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential Roots Upper:', 'blue', attrs=['bold']),potentialRoot_Upper_Sum.getInfo())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "presentTGB_Lower_Abs = presentAGB_Lower_Abs.multiply(rootShootRatioLower.add(1))\n",
    "presentTGB_Upper_Abs = presentAGB_Upper_Abs.multiply(rootShootRatioUpper.add(1))\n",
    "\n",
    "potentialTGB_Lower_Abs = potentialAGB_Lower_Abs.multiply(rootShootRatioLower.add(1))\n",
    "potentialTGB_Upper_Abs = potentialAGB_Upper_Abs.multiply(rootShootRatioUpper.add(1))\n",
    "\n",
    "# presentTGB_Lower_Sum = presentTGB_Lower_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present TGB Lower:', 'blue', attrs=['bold']),presentTGB_Lower_Sum.getInfo())\n",
    "\n",
    "# presentTGB_Upper_Sum = presentTGB_Upper_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present TGB Upper:', 'blue', attrs=['bold']),presentTGB_Upper_Sum.getInfo())\n",
    "\n",
    "# potentialTGB_Lower_Sum = potentialTGB_Lower_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential TGB Lower:', 'blue', attrs=['bold']),potentialTGB_Lower_Sum.getInfo())\n",
    "# potentialTGB_Upper_Sum = potentialTGB_Upper_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential TGB Upper:', 'blue', attrs=['bold']),potentialTGB_Upper_Sum.getInfo())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5 Calculate the Upper and Lower of present and potential Root and PGB"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# load the dead wood and litter layer\n",
    "deadWoodLitterRatio = ee.Image(\"users/leonidmoore/ForestBiomass/DeadWoodLitter/DeadWood_Litter_Ratio_Map\").unmask()\n",
    "deadWoodLitterRatioLower = ee.Image(\"users/leonidmoore/ForestBiomass/DeadWoodLitter/DeadWood_Litter_Ratio_Lower_Map\").unmask()\n",
    "deadWoodLitterRatioUpper = ee.Image(\"users/leonidmoore/ForestBiomass/DeadWoodLitter/DeadWood_Litter_Ratio_Upper_Map\").unmask()\n",
    "# calculate the present and potential PGB\n",
    "presentPGB_Abs = presentTGB_Abs.multiply(deadWoodLitterRatio)\n",
    "potentialPGB_Abs = potentialTGB_Abs.multiply(deadWoodLitterRatio)\n",
    "\n",
    "# calculate the present and potential dead wood and litter\n",
    "presentLitter_Abs = presentTGB_Abs.multiply(deadWoodLitterRatio.subtract(1))\n",
    "potentialLitter_Abs = potentialTGB_Abs.multiply(deadWoodLitterRatio.subtract(1))\n",
    "\n",
    "# calculate the present Dead wood and litter\n",
    "presentLitter_Lower_Abs = presentTGB_Lower_Abs.multiply(deadWoodLitterRatioLower.subtract(1))\n",
    "presentLitter_Upper_Abs = presentTGB_Upper_Abs.multiply(deadWoodLitterRatioUpper.subtract(1))\n",
    "# calculate the potential dead wood and litter\n",
    "potentialLitter_Lower_Abs = potentialTGB_Lower_Abs.multiply(deadWoodLitterRatioLower.subtract(1))\n",
    "potentialLitter_Upper_Abs = potentialTGB_Upper_Abs.multiply(deadWoodLitterRatioUpper.subtract(1))\n",
    "# get the densities\n",
    "presentPGB_Lower = presentTGB_Lower.multiply(deadWoodLitterRatioLower)\n",
    "presentPGB_Upper = presentTGB_Upper.multiply(deadWoodLitterRatioUpper)\n",
    "\n",
    "potentialPGB_Lower = potentialTGB_Lower.multiply(deadWoodLitterRatioLower)\n",
    "potentialPGB_Lower = potentialTGB_Lower.multiply(deadWoodLitterRatioUpper)\n",
    "\n",
    "presentLitter_Lower = presentTGB_Lower.multiply(deadWoodLitterRatioLower.subtract(1))\n",
    "presentLitter_Upper = presentTGB_Upper.multiply(deadWoodLitterRatioUpper.subtract(1))\n",
    "\n",
    "potentialLitter_Lower = potentialTGB_Lower.multiply(deadWoodLitterRatioLower.subtract(1))\n",
    "potentialLitter_Lower = potentialTGB_Lower.multiply(deadWoodLitterRatioUpper.subtract(1))\n",
    "\n",
    "# presentPGB_Abs_Sum = presentPGB_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present PGB:', 'blue', attrs=['bold']),presentPGB_Abs_Sum.getInfo())\n",
    "\n",
    "# potentialPGB_Abs_Sum = potentialPGB_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential PGB:', 'blue', attrs=['bold']),potentialPGB_Abs_Sum.getInfo())\n",
    "\n",
    "# presentLitter_Abs_Sum = presentLitter_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present Dead wood and litter:', 'blue', attrs=['bold']),presentLitter_Abs_Sum.getInfo())\n",
    "\n",
    "# potentialLitter_Abs_Sum = potentialLitter_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential Dead wood and litter:', 'blue', attrs=['bold']),potentialLitter_Abs_Sum.getInfo())\n",
    "\n",
    "# presentLitter_Lower_Sum = presentLitter_Lower_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present Dead wood and litter Lower:', 'blue', attrs=['bold']),presentLitter_Lower_Sum.getInfo())\n",
    "\n",
    "# presentLitter_Upper_Sum = presentLitter_Upper_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present Dead wood and litter Upper:', 'blue', attrs=['bold']),presentLitter_Upper_Sum.getInfo())\n",
    "\n",
    "# potentialLitter_Lower_Sum = potentialLitter_Lower_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential Dead wood and litter Lower:', 'blue', attrs=['bold']),potentialLitter_Lower_Sum.getInfo())\n",
    "\n",
    "# potentialLitter_Upper_Sum = potentialLitter_Upper_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "# # print the estimation out\n",
    "# print(colored('The potential Dead wood and litter Upper:', 'blue', attrs=['bold']),potentialLitter_Upper_Sum.getInfo())\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "# calculate the present PGB Lower and Upper\n",
    "presentPGB_Lower_Abs = presentTGB_Lower_Abs.multiply(deadWoodLitterRatioLower)\n",
    "presentPGB_Upper_Abs = presentTGB_Upper_Abs.multiply(deadWoodLitterRatioUpper)\n",
    "# calculate the potential PGB Lower and Upper\n",
    "potentialPGB_Lower_Abs = potentialTGB_Lower_Abs.multiply(deadWoodLitterRatioLower)\n",
    "potentialPGB_Upper_Abs = potentialTGB_Upper_Abs.multiply(deadWoodLitterRatioUpper)\n",
    "\n",
    "presentPGB_D = presentAGB_Density.multiply(rootShootRatio.add(1)).multiply(deadWoodLitterRatio)\n",
    "potentialPGB_D = potentialAGB_Density.multiply(rootShootRatio.add(1)).multiply(deadWoodLitterRatio)\n",
    "\n",
    "# presentPGB_Lower_Sum = presentPGB_Lower_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present PGB Lower:', 'blue', attrs=['bold']),presentPGB_Lower_Sum.getInfo())\n",
    "\n",
    "# presentPGB_Upper_Sum = presentPGB_Upper_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The present PGB Upper:', 'blue', attrs=['bold']),presentPGB_Upper_Sum.getInfo())\n",
    "\n",
    "# potentialPGB_Lower_Sum = potentialPGB_Lower_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential PGB Lower:', 'blue', attrs=['bold']),potentialPGB_Lower_Sum.getInfo())\n",
    "\n",
    "# potentialPGB_Upper_Sum = potentialPGB_Upper_Abs.multiply(biomeMask).reduceRegion(reducer = ee.Reducer.sum(),\n",
    "#                                                        geometry = unboundedGeo,\n",
    "#                                                        crs = 'EPSG:4326',\n",
    "#                                                        crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "#                                                        maxPixels = 1e9)\n",
    "\n",
    "# # print the estimation out\n",
    "# print(colored('The potential PGB Upper:', 'blue', attrs=['bold']),potentialPGB_Upper_Sum.getInfo())\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6 Export the upper and lower images to Assets"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "# load the carbon density layers\n",
    "SandermannCarbonDiff = ee.Image(\"users/leonidmoore/ForestBiomass/SoilOrganicCarbonModel/SOCS_0_200cm_Diff_1km_Present_subtract_NoLU\").unmask()\n",
    "SandermannCarbonPresent = ee.Image(\"users/leonidmoore/ForestBiomass/SoilOrganicCarbonModel/SOCS_0_200cm_1km_Present\").unmask()\n",
    "\n",
    "# mask the diffrence layer\n",
    "SandermannCarbonLoss = SandermannCarbonDiff.multiply(SandermannCarbonDiff.gt(0))\n",
    "\n",
    "# load the present and potential forest cover\n",
    "presentForestCover = compositeImage.select('PresentTreeCover').unmask() # uniform with potential in the  0-1 scale\n",
    "potentialCoverAdjusted = ee.Image(\"users/leonidmoore/ForestBiomass/Bastin_et_al_2019_Potential_Forest_Cover_Adjusted\").unmask().rename('PotentialForestCover')\n",
    "# define the present and potential forest cover masks\n",
    "presentMask = presentForestCover.gt(0)\n",
    "potentialMask = potentialCoverAdjusted.gte(0.1)\n",
    "\n",
    "# calculate the sum of the potential in soil with the consideration of forest cover\n",
    "SandermannCarbonStockLoss = SandermannCarbonLoss.multiply(pixelArea).divide(1000000000).mask(biomeMask).mask(potentialMask).multiply(potentialCoverAdjusted)\n",
    "\n",
    "# add the soil into the PGB as the total potential\n",
    "potentialTotal_Abs = potentialPGB_Abs.add(SandermannCarbonStockLoss)\n",
    "# compose those bands into an image\n",
    "lowerUpperImage = presentAGB_Lower_Abs.rename('preAGB_Lower').addBands(presentAGB_Upper_Abs.rename('preAGB_Upper')).addBands(potentialAGB_Lower_Abs.rename('potAGB_Lower')).addBands(potentialAGB_Upper_Abs.rename('potAGB_Upper')).addBands(presentRoot_Lower_Abs.rename('preRoot_Lower')).addBands(presentRoot_Upper_Abs.rename('preRoot_Upper')).addBands(potentialRoot_Lower_Abs.rename('potRoot_Lower')).addBands(potentialRoot_Upper_Abs.rename('potRoot_Upper')).addBands(presentLitter_Lower_Abs.rename('preLitter_Lower')).addBands(presentLitter_Upper_Abs.rename('preLitter_Upper')).addBands(potentialLitter_Lower_Abs.rename('potLitter_Lower')).addBands(potentialLitter_Upper_Abs.rename('potLitter_Upper')).addBands(potentialTotal_Abs.rename('PotentialTotal'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'state': 'READY',\n",
       " 'description': 'WK2_Upper_Lower_Uncertainty_maps_Export',\n",
       " 'creation_timestamp_ms': 1690807865941,\n",
       " 'update_timestamp_ms': 1690807865941,\n",
       " 'start_timestamp_ms': 0,\n",
       " 'task_type': 'EXPORT_IMAGE',\n",
       " 'id': 'G3C3PT4WDMMTNANBV4SZ42EC',\n",
       " 'name': 'projects/earthengine-legacy/operations/G3C3PT4WDMMTNANBV4SZ42EC'}"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "exportUpperLower = ee.batch.Export.image.toAsset(image = lowerUpperImage,\n",
    "                                               description = 'WK2_Upper_Lower_Uncertainty_maps_Export',\n",
    "                                               assetId = 'users/leonidmoore/ForestBiomass/UncertaintyFigure/WK2_Lower_Upper_Map',\n",
    "                                               region = unboundedGeo,\n",
    "                                               crs = 'EPSG:4326',\n",
    "                                               crsTransform = [0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "                                               maxPixels = 1e13)\n",
    "# start the export task\n",
    "exportUpperLower.start()\n",
    "# show the task status\n",
    "exportUpperLower.status()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 7 Calculate the Abs for different parts at biome level"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[1m\u001b[34mThe biomass partition results in biome: \n",
      "\u001b[0m\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>PresentAGB</th>\n",
       "      <th>PotentialAGB</th>\n",
       "      <th>PresentRoot</th>\n",
       "      <th>PotentialRoot</th>\n",
       "      <th>PresentTGB</th>\n",
       "      <th>PotentialTGB</th>\n",
       "      <th>PresentLitter</th>\n",
       "      <th>PotentialLitter</th>\n",
       "      <th>PresentPGB</th>\n",
       "      <th>PotentialPGB</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>180.7</td>\n",
       "      <td>227.5</td>\n",
       "      <td>48.0</td>\n",
       "      <td>59.6</td>\n",
       "      <td>228.7</td>\n",
       "      <td>287.2</td>\n",
       "      <td>50.3</td>\n",
       "      <td>63.2</td>\n",
       "      <td>279.0</td>\n",
       "      <td>350.4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>7.9</td>\n",
       "      <td>18.3</td>\n",
       "      <td>2.0</td>\n",
       "      <td>4.7</td>\n",
       "      <td>9.9</td>\n",
       "      <td>23.0</td>\n",
       "      <td>2.2</td>\n",
       "      <td>5.1</td>\n",
       "      <td>12.0</td>\n",
       "      <td>28.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2.8</td>\n",
       "      <td>5.2</td>\n",
       "      <td>0.7</td>\n",
       "      <td>1.3</td>\n",
       "      <td>3.5</td>\n",
       "      <td>6.6</td>\n",
       "      <td>0.8</td>\n",
       "      <td>1.4</td>\n",
       "      <td>4.2</td>\n",
       "      <td>8.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>39.2</td>\n",
       "      <td>61.4</td>\n",
       "      <td>10.0</td>\n",
       "      <td>15.5</td>\n",
       "      <td>49.2</td>\n",
       "      <td>77.0</td>\n",
       "      <td>16.2</td>\n",
       "      <td>25.4</td>\n",
       "      <td>65.5</td>\n",
       "      <td>102.4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>19.3</td>\n",
       "      <td>23.4</td>\n",
       "      <td>5.1</td>\n",
       "      <td>6.2</td>\n",
       "      <td>24.4</td>\n",
       "      <td>29.6</td>\n",
       "      <td>8.0</td>\n",
       "      <td>9.8</td>\n",
       "      <td>32.5</td>\n",
       "      <td>39.4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>35.7</td>\n",
       "      <td>44.7</td>\n",
       "      <td>9.6</td>\n",
       "      <td>12.0</td>\n",
       "      <td>45.3</td>\n",
       "      <td>56.7</td>\n",
       "      <td>36.2</td>\n",
       "      <td>45.3</td>\n",
       "      <td>81.5</td>\n",
       "      <td>102.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>44.3</td>\n",
       "      <td>93.7</td>\n",
       "      <td>13.5</td>\n",
       "      <td>28.4</td>\n",
       "      <td>57.9</td>\n",
       "      <td>122.1</td>\n",
       "      <td>12.7</td>\n",
       "      <td>26.9</td>\n",
       "      <td>70.6</td>\n",
       "      <td>148.9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>3.7</td>\n",
       "      <td>22.5</td>\n",
       "      <td>1.1</td>\n",
       "      <td>6.6</td>\n",
       "      <td>4.8</td>\n",
       "      <td>29.1</td>\n",
       "      <td>1.6</td>\n",
       "      <td>9.6</td>\n",
       "      <td>6.4</td>\n",
       "      <td>38.7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>2.2</td>\n",
       "      <td>3.8</td>\n",
       "      <td>0.7</td>\n",
       "      <td>1.1</td>\n",
       "      <td>2.8</td>\n",
       "      <td>5.0</td>\n",
       "      <td>0.6</td>\n",
       "      <td>1.1</td>\n",
       "      <td>3.4</td>\n",
       "      <td>6.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>3.4</td>\n",
       "      <td>10.1</td>\n",
       "      <td>1.0</td>\n",
       "      <td>2.7</td>\n",
       "      <td>4.4</td>\n",
       "      <td>12.8</td>\n",
       "      <td>1.4</td>\n",
       "      <td>4.2</td>\n",
       "      <td>5.8</td>\n",
       "      <td>17.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>3.4</td>\n",
       "      <td>9.0</td>\n",
       "      <td>1.1</td>\n",
       "      <td>2.7</td>\n",
       "      <td>4.5</td>\n",
       "      <td>11.7</td>\n",
       "      <td>3.6</td>\n",
       "      <td>9.3</td>\n",
       "      <td>8.1</td>\n",
       "      <td>21.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>3.7</td>\n",
       "      <td>9.0</td>\n",
       "      <td>1.1</td>\n",
       "      <td>2.8</td>\n",
       "      <td>4.8</td>\n",
       "      <td>11.7</td>\n",
       "      <td>1.0</td>\n",
       "      <td>2.5</td>\n",
       "      <td>5.8</td>\n",
       "      <td>14.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>2.0</td>\n",
       "      <td>29.0</td>\n",
       "      <td>0.6</td>\n",
       "      <td>7.3</td>\n",
       "      <td>2.7</td>\n",
       "      <td>36.3</td>\n",
       "      <td>0.6</td>\n",
       "      <td>7.6</td>\n",
       "      <td>3.2</td>\n",
       "      <td>43.9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13</th>\n",
       "      <td>1.8</td>\n",
       "      <td>2.7</td>\n",
       "      <td>0.4</td>\n",
       "      <td>0.6</td>\n",
       "      <td>2.2</td>\n",
       "      <td>3.2</td>\n",
       "      <td>0.5</td>\n",
       "      <td>0.7</td>\n",
       "      <td>2.6</td>\n",
       "      <td>4.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>sum</th>\n",
       "      <td>350.1</td>\n",
       "      <td>560.3</td>\n",
       "      <td>94.9</td>\n",
       "      <td>151.5</td>\n",
       "      <td>445.1</td>\n",
       "      <td>712.0</td>\n",
       "      <td>135.7</td>\n",
       "      <td>212.1</td>\n",
       "      <td>580.6</td>\n",
       "      <td>924.2</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "     PresentAGB  PotentialAGB  PresentRoot  PotentialRoot  PresentTGB  \\\n",
       "0         180.7         227.5         48.0           59.6       228.7   \n",
       "1           7.9          18.3          2.0            4.7         9.9   \n",
       "2           2.8           5.2          0.7            1.3         3.5   \n",
       "3          39.2          61.4         10.0           15.5        49.2   \n",
       "4          19.3          23.4          5.1            6.2        24.4   \n",
       "5          35.7          44.7          9.6           12.0        45.3   \n",
       "6          44.3          93.7         13.5           28.4        57.9   \n",
       "7           3.7          22.5          1.1            6.6         4.8   \n",
       "8           2.2           3.8          0.7            1.1         2.8   \n",
       "9           3.4          10.1          1.0            2.7         4.4   \n",
       "10          3.4           9.0          1.1            2.7         4.5   \n",
       "11          3.7           9.0          1.1            2.8         4.8   \n",
       "12          2.0          29.0          0.6            7.3         2.7   \n",
       "13          1.8           2.7          0.4            0.6         2.2   \n",
       "sum       350.1         560.3         94.9          151.5       445.1   \n",
       "\n",
       "     PotentialTGB  PresentLitter  PotentialLitter  PresentPGB  PotentialPGB  \n",
       "0           287.2           50.3             63.2       279.0         350.4  \n",
       "1            23.0            2.2              5.1        12.0          28.1  \n",
       "2             6.6            0.8              1.4         4.2           8.0  \n",
       "3            77.0           16.2             25.4        65.5         102.4  \n",
       "4            29.6            8.0              9.8        32.5          39.4  \n",
       "5            56.7           36.2             45.3        81.5         102.1  \n",
       "6           122.1           12.7             26.9        70.6         148.9  \n",
       "7            29.1            1.6              9.6         6.4          38.7  \n",
       "8             5.0            0.6              1.1         3.4           6.1  \n",
       "9            12.8            1.4              4.2         5.8          17.0  \n",
       "10           11.7            3.6              9.3         8.1          21.0  \n",
       "11           11.7            1.0              2.5         5.8          14.2  \n",
       "12           36.3            0.6              7.6         3.2          43.9  \n",
       "13            3.2            0.5              0.7         2.6           4.0  \n",
       "sum         712.0          135.7            212.1       580.6         924.2  "
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Stack the absolute biomass layers into an Image.\n",
    "absImage = presentAGB_Abs.rename('PresentAGB').addBands(potentialAGB_Abs.rename('PotentialAGB')).addBands(presentRoot_Abs.rename('PresentRoot')).addBands(potentialRoot_Abs.rename('PotentialRoot')).addBands(presentTGB_Abs.rename('PresentTGB')).addBands(potentialTGB_Abs.rename('PotentialTGB')).addBands(presentLitter_Abs.rename('PresentLitter')).addBands(potentialLitter_Abs.rename('PotentialLitter')).addBands(presentPGB_Abs.rename('PresentPGB')).addBands(potentialPGB_Abs.rename('PotentialPGB'))\n",
    "\n",
    "# define the function to do the biome level statistics which could be applied by map      \n",
    "def biomeLevelStat(biome):\n",
    "    perBiomeMask = biomeLayer.eq(ee.Number(biome))\n",
    "    masked_img = absImage.mask(perBiomeMask)\n",
    "    output = masked_img.reduceRegion(reducer= ee.Reducer.sum(),\n",
    "                                     geometry= unboundedGeo,\n",
    "                                     crs='EPSG:4326',\n",
    "                                     crsTransform=[0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "                                     maxPixels= 1e13)\n",
    "    return output#.getInfo().get('Present')\n",
    "\n",
    "\n",
    "biomeList = ee.List([1,2,3,4,5,6,7,8,9,10,11,12,13,14])\n",
    "statisticTable = biomeList.map(biomeLevelStat).getInfo()\n",
    "# transform into data frame\n",
    "outputTable = pd.DataFrame(statisticTable,columns =['PresentAGB','PotentialAGB','PresentRoot','PotentialRoot','PresentTGB','PotentialTGB','PresentLitter','PotentialLitter','PresentPGB','PotentialPGB']).round(1)\n",
    "outputTable.loc['sum'] = outputTable.sum() \n",
    "outputTable.to_csv('Data/BiomeLevelStatistics/StatisticsForModels/WK2_Abs_for_diff_parts_at_Biome_Level.csv',header=True,mode='w+')\n",
    "print(colored('The biomass partition results in biome: \\n', 'blue', attrs=['bold']))\n",
    "outputTable.head(15)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [],
   "source": [
    "# If you got the error 'EEException: Too many concurrent aggregations.', please re-run this chunck of code again."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 7 Calculate the Abs for different parts at biome level"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[1m\u001b[34mThe biomass partition results in biome: \n",
      "\u001b[0m\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>preAGB_Lower</th>\n",
       "      <th>preAGB_Upper</th>\n",
       "      <th>potAGB_Lower</th>\n",
       "      <th>potAGB_Upper</th>\n",
       "      <th>preRoot_Lower</th>\n",
       "      <th>preRoot_Upper</th>\n",
       "      <th>potRoot_Lower</th>\n",
       "      <th>potRoot_Upper</th>\n",
       "      <th>preLitter_Lower</th>\n",
       "      <th>preLitter_Upper</th>\n",
       "      <th>potLitter_Lower</th>\n",
       "      <th>potLitter_Upper</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>179.4</td>\n",
       "      <td>181.9</td>\n",
       "      <td>204.7</td>\n",
       "      <td>235.6</td>\n",
       "      <td>37.6</td>\n",
       "      <td>59.1</td>\n",
       "      <td>42.7</td>\n",
       "      <td>74.9</td>\n",
       "      <td>32.6</td>\n",
       "      <td>72.3</td>\n",
       "      <td>37.1</td>\n",
       "      <td>93.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>7.8</td>\n",
       "      <td>7.9</td>\n",
       "      <td>15.6</td>\n",
       "      <td>20.9</td>\n",
       "      <td>1.6</td>\n",
       "      <td>2.4</td>\n",
       "      <td>3.2</td>\n",
       "      <td>6.2</td>\n",
       "      <td>1.4</td>\n",
       "      <td>3.1</td>\n",
       "      <td>2.8</td>\n",
       "      <td>8.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2.7</td>\n",
       "      <td>2.8</td>\n",
       "      <td>4.6</td>\n",
       "      <td>5.9</td>\n",
       "      <td>0.6</td>\n",
       "      <td>0.8</td>\n",
       "      <td>1.0</td>\n",
       "      <td>1.7</td>\n",
       "      <td>0.5</td>\n",
       "      <td>1.1</td>\n",
       "      <td>0.8</td>\n",
       "      <td>2.3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>38.8</td>\n",
       "      <td>39.7</td>\n",
       "      <td>51.8</td>\n",
       "      <td>64.7</td>\n",
       "      <td>8.8</td>\n",
       "      <td>11.4</td>\n",
       "      <td>11.7</td>\n",
       "      <td>18.6</td>\n",
       "      <td>14.3</td>\n",
       "      <td>18.9</td>\n",
       "      <td>19.0</td>\n",
       "      <td>30.8</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>19.0</td>\n",
       "      <td>19.5</td>\n",
       "      <td>17.1</td>\n",
       "      <td>22.7</td>\n",
       "      <td>4.4</td>\n",
       "      <td>5.9</td>\n",
       "      <td>4.0</td>\n",
       "      <td>6.9</td>\n",
       "      <td>7.0</td>\n",
       "      <td>9.4</td>\n",
       "      <td>6.3</td>\n",
       "      <td>10.9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>35.4</td>\n",
       "      <td>36.0</td>\n",
       "      <td>37.9</td>\n",
       "      <td>46.4</td>\n",
       "      <td>8.7</td>\n",
       "      <td>10.7</td>\n",
       "      <td>9.4</td>\n",
       "      <td>13.9</td>\n",
       "      <td>30.0</td>\n",
       "      <td>43.9</td>\n",
       "      <td>32.1</td>\n",
       "      <td>56.7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>43.9</td>\n",
       "      <td>44.8</td>\n",
       "      <td>78.1</td>\n",
       "      <td>109.2</td>\n",
       "      <td>11.4</td>\n",
       "      <td>15.7</td>\n",
       "      <td>20.0</td>\n",
       "      <td>38.5</td>\n",
       "      <td>8.3</td>\n",
       "      <td>18.2</td>\n",
       "      <td>14.7</td>\n",
       "      <td>44.3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>3.7</td>\n",
       "      <td>3.8</td>\n",
       "      <td>18.7</td>\n",
       "      <td>26.0</td>\n",
       "      <td>0.9</td>\n",
       "      <td>1.2</td>\n",
       "      <td>4.8</td>\n",
       "      <td>8.9</td>\n",
       "      <td>1.4</td>\n",
       "      <td>1.9</td>\n",
       "      <td>7.0</td>\n",
       "      <td>12.9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>2.1</td>\n",
       "      <td>2.2</td>\n",
       "      <td>3.2</td>\n",
       "      <td>4.4</td>\n",
       "      <td>0.6</td>\n",
       "      <td>0.8</td>\n",
       "      <td>0.8</td>\n",
       "      <td>1.5</td>\n",
       "      <td>0.4</td>\n",
       "      <td>0.9</td>\n",
       "      <td>0.6</td>\n",
       "      <td>1.8</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>3.4</td>\n",
       "      <td>3.5</td>\n",
       "      <td>8.0</td>\n",
       "      <td>12.0</td>\n",
       "      <td>0.8</td>\n",
       "      <td>1.1</td>\n",
       "      <td>1.8</td>\n",
       "      <td>3.7</td>\n",
       "      <td>1.2</td>\n",
       "      <td>1.7</td>\n",
       "      <td>2.9</td>\n",
       "      <td>5.8</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>3.3</td>\n",
       "      <td>3.5</td>\n",
       "      <td>6.9</td>\n",
       "      <td>10.6</td>\n",
       "      <td>0.9</td>\n",
       "      <td>1.3</td>\n",
       "      <td>1.8</td>\n",
       "      <td>3.9</td>\n",
       "      <td>2.9</td>\n",
       "      <td>4.5</td>\n",
       "      <td>5.9</td>\n",
       "      <td>13.6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>3.7</td>\n",
       "      <td>3.8</td>\n",
       "      <td>7.0</td>\n",
       "      <td>10.4</td>\n",
       "      <td>0.9</td>\n",
       "      <td>1.5</td>\n",
       "      <td>1.8</td>\n",
       "      <td>4.2</td>\n",
       "      <td>0.1</td>\n",
       "      <td>2.1</td>\n",
       "      <td>0.2</td>\n",
       "      <td>5.8</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>2.0</td>\n",
       "      <td>2.1</td>\n",
       "      <td>23.0</td>\n",
       "      <td>36.1</td>\n",
       "      <td>0.5</td>\n",
       "      <td>0.8</td>\n",
       "      <td>5.0</td>\n",
       "      <td>11.2</td>\n",
       "      <td>0.1</td>\n",
       "      <td>1.1</td>\n",
       "      <td>0.6</td>\n",
       "      <td>18.9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13</th>\n",
       "      <td>1.8</td>\n",
       "      <td>1.8</td>\n",
       "      <td>2.3</td>\n",
       "      <td>2.9</td>\n",
       "      <td>0.3</td>\n",
       "      <td>0.5</td>\n",
       "      <td>0.4</td>\n",
       "      <td>0.9</td>\n",
       "      <td>0.3</td>\n",
       "      <td>0.7</td>\n",
       "      <td>0.4</td>\n",
       "      <td>1.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>sum</th>\n",
       "      <td>347.0</td>\n",
       "      <td>353.3</td>\n",
       "      <td>478.9</td>\n",
       "      <td>607.8</td>\n",
       "      <td>78.0</td>\n",
       "      <td>113.2</td>\n",
       "      <td>108.4</td>\n",
       "      <td>195.0</td>\n",
       "      <td>100.5</td>\n",
       "      <td>179.8</td>\n",
       "      <td>130.4</td>\n",
       "      <td>306.2</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "     preAGB_Lower  preAGB_Upper  potAGB_Lower  potAGB_Upper  preRoot_Lower  \\\n",
       "0           179.4         181.9         204.7         235.6           37.6   \n",
       "1             7.8           7.9          15.6          20.9            1.6   \n",
       "2             2.7           2.8           4.6           5.9            0.6   \n",
       "3            38.8          39.7          51.8          64.7            8.8   \n",
       "4            19.0          19.5          17.1          22.7            4.4   \n",
       "5            35.4          36.0          37.9          46.4            8.7   \n",
       "6            43.9          44.8          78.1         109.2           11.4   \n",
       "7             3.7           3.8          18.7          26.0            0.9   \n",
       "8             2.1           2.2           3.2           4.4            0.6   \n",
       "9             3.4           3.5           8.0          12.0            0.8   \n",
       "10            3.3           3.5           6.9          10.6            0.9   \n",
       "11            3.7           3.8           7.0          10.4            0.9   \n",
       "12            2.0           2.1          23.0          36.1            0.5   \n",
       "13            1.8           1.8           2.3           2.9            0.3   \n",
       "sum         347.0         353.3         478.9         607.8           78.0   \n",
       "\n",
       "     preRoot_Upper  potRoot_Lower  potRoot_Upper  preLitter_Lower  \\\n",
       "0             59.1           42.7           74.9             32.6   \n",
       "1              2.4            3.2            6.2              1.4   \n",
       "2              0.8            1.0            1.7              0.5   \n",
       "3             11.4           11.7           18.6             14.3   \n",
       "4              5.9            4.0            6.9              7.0   \n",
       "5             10.7            9.4           13.9             30.0   \n",
       "6             15.7           20.0           38.5              8.3   \n",
       "7              1.2            4.8            8.9              1.4   \n",
       "8              0.8            0.8            1.5              0.4   \n",
       "9              1.1            1.8            3.7              1.2   \n",
       "10             1.3            1.8            3.9              2.9   \n",
       "11             1.5            1.8            4.2              0.1   \n",
       "12             0.8            5.0           11.2              0.1   \n",
       "13             0.5            0.4            0.9              0.3   \n",
       "sum          113.2          108.4          195.0            100.5   \n",
       "\n",
       "     preLitter_Upper  potLitter_Lower  potLitter_Upper  \n",
       "0               72.3             37.1             93.2  \n",
       "1                3.1              2.8              8.1  \n",
       "2                1.1              0.8              2.3  \n",
       "3               18.9             19.0             30.8  \n",
       "4                9.4              6.3             10.9  \n",
       "5               43.9             32.1             56.7  \n",
       "6               18.2             14.7             44.3  \n",
       "7                1.9              7.0             12.9  \n",
       "8                0.9              0.6              1.8  \n",
       "9                1.7              2.9              5.8  \n",
       "10               4.5              5.9             13.6  \n",
       "11               2.1              0.2              5.8  \n",
       "12               1.1              0.6             18.9  \n",
       "13               0.7              0.4              1.1  \n",
       "sum            179.8            130.4            306.2  "
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Stack the absolute biomass layers into an Image.\n",
    "absPotentialImage = presentAGB_Lower_Abs.rename('preAGB_Lower').addBands(presentAGB_Upper_Abs.rename('preAGB_Upper')).addBands(potentialAGB_Lower_Abs.rename('potAGB_Lower')).addBands(potentialAGB_Upper_Abs.rename('potAGB_Upper')).addBands(presentRoot_Lower_Abs.rename('preRoot_Lower')).addBands(presentRoot_Upper_Abs.rename('preRoot_Upper')).addBands(potentialRoot_Lower_Abs.rename('potRoot_Lower')).addBands(potentialRoot_Upper_Abs.rename('potRoot_Upper')).addBands(presentLitter_Lower_Abs.rename('preLitter_Lower')).addBands(presentLitter_Upper_Abs.rename('preLitter_Upper')).addBands(potentialLitter_Lower_Abs.rename('potLitter_Lower')).addBands(potentialLitter_Upper_Abs.rename('potLitter_Upper'))\n",
    "# define the function to do the biome level statistics which could be applied by map      \n",
    "def biomeLevelStat(biome):\n",
    "    perBiomeMask = biomeLayer.eq(ee.Number(biome))\n",
    "    masked_img = absPotentialImage.mask(perBiomeMask)\n",
    "    output = masked_img.reduceRegion(reducer= ee.Reducer.sum(),\n",
    "                                     geometry= unboundedGeo,\n",
    "                                     crs='EPSG:4326',\n",
    "                                     crsTransform=[0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "                                     maxPixels= 1e13)\n",
    "    return output#.getInfo().get('Present')\n",
    "\n",
    "\n",
    "biomeList = ee.List([1,2,3,4,5,6,7,8,9,10,11,12,13,14])\n",
    "statisticTable = biomeList.map(biomeLevelStat).getInfo()\n",
    "# transform into data frame\n",
    "outputTable = pd.DataFrame(statisticTable,columns =['preAGB_Lower','preAGB_Upper','potAGB_Lower','potAGB_Upper','preRoot_Lower','preRoot_Upper','potRoot_Lower','potRoot_Upper','preLitter_Lower','preLitter_Upper','potLitter_Lower','potLitter_Upper']).round(1)\n",
    "outputTable.loc['sum'] = outputTable.sum() \n",
    "outputTable.to_csv('Data/BiomeLevelStatistics/StatisticsForModels/WK2_Uncertainty_for_diff_parts_at_Biome_Level.csv',header=True,mode='w+')\n",
    "print(colored('The biomass partition results in biome: \\n', 'blue', attrs=['bold']))\n",
    "outputTable.head(15)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# If you got the error 'EEException: Too many concurrent aggregations.', please re-run this chunck of code again."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 8 Calculate the uncertainty for each land use types"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Load all the landuse type layers\n",
    "croplandOrg = ee.Image(\"users/leonidmoore/ForestBiomass/HYDE31/cropland_Percent\").rename('cropland').divide(100).reproject(crs=stdProj);\n",
    "grazingOrg = ee.Image(\"users/leonidmoore/ForestBiomass/HYDE31/grazing_Percent\").rename('grazing').divide(100).reproject(crs=stdProj);\n",
    "pastureOrg = ee.Image(\"users/leonidmoore/ForestBiomass/HYDE31/pasture_Percent\").rename('pasture').divide(100).reproject(crs=stdProj);\n",
    "rangelandOrg = ee.Image(\"users/leonidmoore/ForestBiomass/HYDE31/rangeland_Percent\").rename('rangeland').divide(100).reproject(crs=stdProj);\n",
    "urbanOrg = compositeImage.select(['LandCoverClass_Urban_Builtup']).divide(100).unmask().reproject(crs=stdProj);\n",
    "snowIceOrg = compositeImageNew.select(['ConsensusLandCoverClass_Snow_Ice']).divide(100).unmask().reproject(crs=stdProj);\n",
    "openWaterOrg = compositeImageNew.select(['ConsensusLandCoverClass_Open_Water']).divide(100).unmask().reproject(crs=stdProj);\n",
    "# define the total landcover types\n",
    "sumCover = presentForestCover.add(pastureOrg).add(rangelandOrg).add(croplandOrg).add(urbanOrg).add(openWaterOrg).add(snowIceOrg);\n",
    "oneSubtract = ee.Image(1).subtract(sumCover);\n",
    "freeland = oneSubtract.multiply(oneSubtract.gte(0));\n",
    "# get the scale ratio for pixels with sumCover larger than 1\n",
    "scaleRatio = ee.Image(1).subtract(presentForestCover).multiply(oneSubtract.lt(0)).divide(sumCover.subtract(presentForestCover).multiply(oneSubtract.lt(0))).unmask();\n",
    "# get the ratio of these three disturbed maps\n",
    "pasture = pastureOrg.multiply(scaleRatio).multiply(oneSubtract.lt(0)).add(pastureOrg.multiply(oneSubtract.gte(0))).unmask();\n",
    "rangeland = rangelandOrg.multiply(scaleRatio).multiply(oneSubtract.lt(0)).add(rangelandOrg.multiply(oneSubtract.gte(0))).unmask();\n",
    "cropland = croplandOrg.multiply(scaleRatio).multiply(oneSubtract.lt(0)).add(croplandOrg.multiply(oneSubtract.gte(0))).unmask();\n",
    "urban = urbanOrg.multiply(scaleRatio).multiply(oneSubtract.lt(0)).add(urbanOrg.multiply(oneSubtract.gte(0))).unmask();\n",
    "openWater = openWaterOrg.multiply(scaleRatio).multiply(oneSubtract.lt(0)).add(openWaterOrg.multiply(oneSubtract.gte(0))).unmask();\n",
    "snowIce = snowIceOrg.multiply(scaleRatio).multiply(oneSubtract.lt(0)).add(snowIceOrg.multiply(oneSubtract.gte(0))).unmask();\n",
    "sumTT = presentForestCover.add(pasture).add(rangeland).add(cropland).add(urban).add(freeland).add(openWater).add(snowIce).unmask();\n",
    "\n",
    "effectivePotentialMask = freeland.add(rangeland).add(pasture).add(cropland).add(urban).gt(0);\n",
    "# there are some pixels without any landcover survived but with open water and ice and snow. here we mask these pixels out\n",
    "sumlandCover = pastureOrg.add(rangelandOrg).add(croplandOrg).add(urbanOrg).add(freeland);\n",
    "restorationMap = potentialForestCover.subtract(presentForestCover).mask(sumlandCover.neq(0)).unmask();\n",
    "\n",
    "# sum all these scaled layersv\n",
    "scaledSum = pasture.add(rangeland).add(cropland).add(urban).add(freeland);\n",
    "potentialCoverFinal = restorationMap.add(presentForestCover);\n",
    "# allocate the potential equally to each layer\n",
    "freelandPotentialCover = freeland.divide(scaledSum).multiply(restorationMap).unmask();\n",
    "rangelandPotentialCover = rangeland.divide(scaledSum).multiply(restorationMap).unmask();\n",
    "pasturePotentialCover = pasture.divide(scaledSum).multiply(restorationMap).unmask();\n",
    "croplandPotentialCover = cropland.divide(scaledSum).multiply(restorationMap).unmask();\n",
    "urbanPotentialCover = urban.divide(scaledSum).multiply(restorationMap).unmask();\n",
    "#  allocate the freeland potential in pixels with forest cover larger than 10% to conservation potential\n",
    "freelandForConsevation = freelandPotentialCover.multiply(presentForestCover.gte(0.1)).unmask();\n",
    "maximumPotentialCover = freelandForConsevation.add(presentForestCover);\n",
    "# calucate the reall freeland outside of forest\n",
    "freelandLeftMap = freelandPotentialCover.subtract(freelandForConsevation)# the left positive pixels are real freeland pixels"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "# define the fucntion that can do the carbon partitionning into differen land use types\n",
    "def absCalculationFunc(presentDensity,potentialDensity):\n",
    "    # calculate the existing carbon, present potential carbon and absolute potential carbon in forests\n",
    "    absoluteImage1 = presentDensity.multiply(pixelArea).multiply(presentMask).divide(1000000000).rename('Present')\n",
    "    absoluteImage3 = potentialDensity.multiply(pixelArea).multiply(potentialMask).divide(1000000000).rename('AbsolutePotential')\n",
    "    # get the sum of the potential covers\n",
    "    potentialCoverSum = freelandLeftMap.add(rangelandPotentialCover).add(pasturePotentialCover).add(croplandPotentialCover).add(urbanPotentialCover)\n",
    "\n",
    "    trueRestorationPotential = absoluteImage3.subtract(absoluteImage1).multiply(1000000000)\n",
    "    ratioPotentialBiomassDensity = absoluteImage1.multiply(potentialForestCover.divide(presentForestCover))\n",
    "    #  get the real for conservation potential\n",
    "    realDensityIncreased = absoluteImage3.subtract(absoluteImage1).mask(absoluteImage3.subtract(ratioPotentialBiomassDensity).gt(0)).unmask()\n",
    "    realDensityNotIncreased = absoluteImage3.subtract(absoluteImage1).mask(absoluteImage3.subtract(ratioPotentialBiomassDensity).lte(0)).unmask()\n",
    "    trueReforestationPotential = realDensityNotIncreased.add(realDensityIncreased.multiply(ee.Image(1).subtract(presentForestCover.divide(potentialForestCover))))\n",
    "\n",
    "    conservationPotentialPart1 = realDensityIncreased.multiply(presentForestCover.add(freelandForConsevation).divide(potentialForestCover))\n",
    "    conservationPotentialPart2 = realDensityNotIncreased.multiply(freelandForConsevation.divide(potentialForestCover.subtract(presentForestCover)))\n",
    "\n",
    "    absoluteImage2 = conservationPotentialPart1.add(conservationPotentialPart2).add(absoluteImage1).rename('PresentPotential')\n",
    "\n",
    "    trueReforestationPotential = absoluteImage3.subtract(absoluteImage2)\n",
    "\n",
    "    absoluteImage4 = trueReforestationPotential.multiply(freelandLeftMap.divide(potentialCoverSum)).rename('FreelandPotential')\n",
    "    absoluteImage5 = trueReforestationPotential.multiply(rangelandPotentialCover.divide(potentialCoverSum)).rename('RangelandPotential')\n",
    "    absoluteImage6 = trueReforestationPotential.multiply(pasturePotentialCover.divide(potentialCoverSum)).rename('PasturePotential')\n",
    "    absoluteImage7 = trueReforestationPotential.multiply(croplandPotentialCover.divide(potentialCoverSum)).rename('CroplandPotential')\n",
    "    absoluteImage8 = trueReforestationPotential.multiply(urbanPotentialCover.divide(potentialCoverSum)).rename('UrbanPotential')\n",
    "    \n",
    "    # Stack the absolute biomass layers into an Image.\n",
    "    absPotentialImage = absoluteImage1.addBands(absoluteImage2).addBands(absoluteImage3).addBands(absoluteImage4).addBands(absoluteImage5).addBands(absoluteImage6).addBands(absoluteImage7).addBands(absoluteImage8)\n",
    "    \n",
    "    output = absPotentialImage.mask(biomeMask).reduceRegion(reducer= ee.Reducer.sum(),\n",
    "                                     geometry= unboundedGeo,\n",
    "                                     crs='EPSG:4326',\n",
    "                                     crsTransform=[0.008333333333333333,0,-180,0,-0.008333333333333333,90],\n",
    "                                     maxPixels= 1e13)\n",
    "    statisticTable = output.getInfo()\n",
    "    # # transform into data frame\n",
    "    outputTable = pd.DataFrame(statisticTable,index=['i',],columns =['Types','Present','PresentPotential','AbsolutePotential','FreelandPotential','RangelandPotential','PasturePotential','CroplandPotential','UrbanPotential']).round(1)\n",
    "    # return the outputTable\n",
    "    return outputTable"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "# for AGBs\n",
    "resultAGB_Lower = absCalculationFunc(presentDensity = presentAGB_Lower,potentialDensity = potentialAGB_AdjLower)\n",
    "resultAGB_Upper = absCalculationFunc(presentDensity = presentAGB_Upper,potentialDensity = potentialAGB_AdjUpper)\n",
    "\n",
    "# For roots\n",
    "presentRoot_Lower = presentAGB_Lower.multiply(rootShootRatioLower).mask(presentMask)\n",
    "potentialRoot_Lower = potentialAGB_AdjLower.multiply(rootShootRatioLower).mask(potentialMask)\n",
    "presentRoot_Upper = presentAGB_Upper.multiply(rootShootRatioUpper).mask(presentMask)\n",
    "potentialRoot_Upper = potentialAGB_AdjUpper.multiply(rootShootRatioUpper).mask(potentialMask)\n",
    "\n",
    "resultRoot_Lower = absCalculationFunc(presentDensity = presentRoot_Lower,potentialDensity = potentialRoot_Lower)\n",
    "resultRoot_Upper = absCalculationFunc(presentDensity = presentRoot_Upper,potentialDensity = potentialRoot_Upper)\n",
    "# for Litter\n",
    "presentLitter_Lower = presentLitter_Lower_Abs.divide(pixelArea).multiply(1000000000).mask(presentMask)\n",
    "potentialLitter_Lower = potentialLitter_Lower_Abs.divide(pixelArea).multiply(1000000000).mask(potentialMask)\n",
    "presentLitter_Upper = presentLitter_Upper_Abs.divide(pixelArea).multiply(1000000000).mask(presentMask)\n",
    "potentialLitter_Upper = potentialLitter_Upper_Abs.divide(pixelArea).multiply(1000000000).mask(potentialMask)\n",
    "\n",
    "resultLitter_Lower = absCalculationFunc(presentDensity = presentLitter_Lower,potentialDensity = potentialLitter_Lower)\n",
    "resultLitter_Upper = absCalculationFunc(presentDensity = presentLitter_Upper,potentialDensity = potentialLitter_Upper)\n",
    "# for TGB\n",
    "presentTGB_Lower = presentTGB_Lower_Abs.divide(pixelArea).multiply(1000000000).mask(presentMask)\n",
    "potentialTGB_Lower = potentialTGB_Lower_Abs.divide(pixelArea).multiply(1000000000).mask(potentialMask)\n",
    "presentTGB_Upper = presentTGB_Upper_Abs.divide(pixelArea).multiply(1000000000).mask(presentMask)\n",
    "potentialTGB_Upper = potentialTGB_Upper_Abs.divide(pixelArea).multiply(1000000000).mask(potentialMask)\n",
    "\n",
    "resultTGB_Lower = absCalculationFunc(presentDensity = presentTGB_Lower,potentialDensity = potentialTGB_Lower)\n",
    "resultTGB_Upper = absCalculationFunc(presentDensity = presentTGB_Upper,potentialDensity = potentialTGB_Upper)\n",
    "# For PGB\n",
    "presentPGB_Lower = presentPGB_Lower_Abs.divide(pixelArea).multiply(1000000000).mask(presentMask)\n",
    "potentialPGB_Lower = potentialPGB_Lower_Abs.divide(pixelArea).multiply(1000000000).mask(potentialMask)\n",
    "presentPGB_Upper = presentPGB_Upper_Abs.divide(pixelArea).multiply(1000000000).mask(presentMask)\n",
    "potentialPGB_Upper = potentialPGB_Upper_Abs.divide(pixelArea).multiply(1000000000).mask(potentialMask)\n",
    "\n",
    "resultPGB_Lower = absCalculationFunc(presentDensity = presentPGB_Lower,potentialDensity = potentialPGB_Lower)\n",
    "resultPGB_Upper = absCalculationFunc(presentDensity = presentPGB_Upper,potentialDensity = potentialPGB_Upper)\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "# concat all the tables into a data frame\n",
    "concatTable = pd.concat([resultAGB_Lower,resultAGB_Upper,resultRoot_Lower,resultRoot_Upper,resultLitter_Lower,resultLitter_Upper,resultTGB_Lower,resultTGB_Upper,resultPGB_Lower,resultPGB_Upper])\n",
    "\n",
    "finalTable = concatTable.assign(Types=['AGB_Lower',\n",
    "                                       'AGB_Upper',\n",
    "                                       'Root_Lower',\n",
    "                                       'Root_Upper',\n",
    "                                       'Litter_Lower',\n",
    "                                       'Litter_Upper',\n",
    "                                       'TGB_Lower',\n",
    "                                       'TGB_Upper',\n",
    "                                       'PGB_Lower',\n",
    "                                       'PGB_Upper'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[1m\u001b[34mThe biomass partition results in biome: \n",
      "\u001b[0m\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Types</th>\n",
       "      <th>Present</th>\n",
       "      <th>PresentPotential</th>\n",
       "      <th>AbsolutePotential</th>\n",
       "      <th>FreelandPotential</th>\n",
       "      <th>RangelandPotential</th>\n",
       "      <th>PasturePotential</th>\n",
       "      <th>CroplandPotential</th>\n",
       "      <th>UrbanPotential</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>i</th>\n",
       "      <td>AGB_Lower</td>\n",
       "      <td>347.0</td>\n",
       "      <td>392.8</td>\n",
       "      <td>507.1</td>\n",
       "      <td>34.1</td>\n",
       "      <td>23.7</td>\n",
       "      <td>22.7</td>\n",
       "      <td>32.4</td>\n",
       "      <td>1.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>i</th>\n",
       "      <td>AGB_Upper</td>\n",
       "      <td>353.3</td>\n",
       "      <td>443.5</td>\n",
       "      <td>619.6</td>\n",
       "      <td>54.1</td>\n",
       "      <td>37.7</td>\n",
       "      <td>34.1</td>\n",
       "      <td>48.1</td>\n",
       "      <td>1.5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>i</th>\n",
       "      <td>Root_Lower</td>\n",
       "      <td>78.0</td>\n",
       "      <td>88.4</td>\n",
       "      <td>114.7</td>\n",
       "      <td>4.1</td>\n",
       "      <td>2.5</td>\n",
       "      <td>4.5</td>\n",
       "      <td>5.4</td>\n",
       "      <td>0.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>i</th>\n",
       "      <td>Root_Upper</td>\n",
       "      <td>113.2</td>\n",
       "      <td>142.0</td>\n",
       "      <td>198.2</td>\n",
       "      <td>8.8</td>\n",
       "      <td>5.5</td>\n",
       "      <td>9.3</td>\n",
       "      <td>10.9</td>\n",
       "      <td>0.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>i</th>\n",
       "      <td>Litter_Lower</td>\n",
       "      <td>100.4</td>\n",
       "      <td>113.8</td>\n",
       "      <td>140.3</td>\n",
       "      <td>5.6</td>\n",
       "      <td>2.0</td>\n",
       "      <td>4.2</td>\n",
       "      <td>5.9</td>\n",
       "      <td>0.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>i</th>\n",
       "      <td>Litter_Upper</td>\n",
       "      <td>179.7</td>\n",
       "      <td>226.9</td>\n",
       "      <td>312.6</td>\n",
       "      <td>15.7</td>\n",
       "      <td>6.7</td>\n",
       "      <td>12.1</td>\n",
       "      <td>15.2</td>\n",
       "      <td>0.5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>i</th>\n",
       "      <td>TGB_Lower</td>\n",
       "      <td>425.0</td>\n",
       "      <td>481.2</td>\n",
       "      <td>621.8</td>\n",
       "      <td>21.1</td>\n",
       "      <td>12.0</td>\n",
       "      <td>22.9</td>\n",
       "      <td>28.8</td>\n",
       "      <td>0.9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>i</th>\n",
       "      <td>TGB_Upper</td>\n",
       "      <td>466.5</td>\n",
       "      <td>585.5</td>\n",
       "      <td>817.8</td>\n",
       "      <td>35.2</td>\n",
       "      <td>20.5</td>\n",
       "      <td>37.0</td>\n",
       "      <td>45.4</td>\n",
       "      <td>1.3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>i</th>\n",
       "      <td>PGB_Lower</td>\n",
       "      <td>525.4</td>\n",
       "      <td>594.9</td>\n",
       "      <td>762.0</td>\n",
       "      <td>26.7</td>\n",
       "      <td>14.0</td>\n",
       "      <td>27.1</td>\n",
       "      <td>34.7</td>\n",
       "      <td>1.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>i</th>\n",
       "      <td>PGB_Upper</td>\n",
       "      <td>646.2</td>\n",
       "      <td>812.4</td>\n",
       "      <td>1130.4</td>\n",
       "      <td>50.9</td>\n",
       "      <td>27.2</td>\n",
       "      <td>49.0</td>\n",
       "      <td>60.6</td>\n",
       "      <td>1.8</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          Types  Present  PresentPotential  AbsolutePotential  \\\n",
       "i     AGB_Lower    347.0             392.8              507.1   \n",
       "i     AGB_Upper    353.3             443.5              619.6   \n",
       "i    Root_Lower     78.0              88.4              114.7   \n",
       "i    Root_Upper    113.2             142.0              198.2   \n",
       "i  Litter_Lower    100.4             113.8              140.3   \n",
       "i  Litter_Upper    179.7             226.9              312.6   \n",
       "i     TGB_Lower    425.0             481.2              621.8   \n",
       "i     TGB_Upper    466.5             585.5              817.8   \n",
       "i     PGB_Lower    525.4             594.9              762.0   \n",
       "i     PGB_Upper    646.2             812.4             1130.4   \n",
       "\n",
       "   FreelandPotential  RangelandPotential  PasturePotential  CroplandPotential  \\\n",
       "i               34.1                23.7              22.7               32.4   \n",
       "i               54.1                37.7              34.1               48.1   \n",
       "i                4.1                 2.5               4.5                5.4   \n",
       "i                8.8                 5.5               9.3               10.9   \n",
       "i                5.6                 2.0               4.2                5.9   \n",
       "i               15.7                 6.7              12.1               15.2   \n",
       "i               21.1                12.0              22.9               28.8   \n",
       "i               35.2                20.5              37.0               45.4   \n",
       "i               26.7                14.0              27.1               34.7   \n",
       "i               50.9                27.2              49.0               60.6   \n",
       "\n",
       "   UrbanPotential  \n",
       "i             1.0  \n",
       "i             1.5  \n",
       "i             0.1  \n",
       "i             0.2  \n",
       "i             0.2  \n",
       "i             0.5  \n",
       "i             0.9  \n",
       "i             1.3  \n",
       "i             1.1  \n",
       "i             1.8  "
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "finalTable.to_csv('BiomeLevelStatistics/StatisticsForModels/WK2_Landuse_type_Uncertainty.csv',header=True,mode='w+')\n",
    "print(colored('The biomass partition results in biome: \\n', 'blue', attrs=['bold']))\n",
    "finalTable.head(10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# If you got the error 'EEException: Too many concurrent aggregations.', please re-run this chunck of code again."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
